Traffic matrix estimation method and apparatus

ABSTRACT

A method and apparatus for the estimation of traffic matrices in a network are disclosed.

This application is a continuation of, and claims priority to, priorapplication Ser. No. 10/627,767 filed Jul. 25, 2003 now U.S. Pat. No.7,293,086, which claims the benefit of Provisional Application60/398,474 filed Jul. 25, 2002, both of which are incorporated herein byreference.

BACKGROUND OF THE INVENTION

The present invention relates generally to network management and, moreparticularly, to estimation of traffic matrices in a network.

A traffic matrix provides, for every ingress point i into a network andegress: point j out of the network, the volume of traffic T_(i,j) from ito j over a given time interval. The traffic matrix, together withnetwork topology and routing and fault data, can facilitate diagnosisand management of network congestion—as well as provide critical inputsto network design, capacity planning and business planning.Unfortunately, traffic matrices are generally unavailable in largeoperational Internet Protocol (“IP”) networks. Rather, typicalproduction systems gather data on resource utilization at network nodesand links (e.g. link loads); end-to-end performance metrics for specifictransactions (such as one way delay statistics for packets exchangedbetween measurement servers at the network edge); and status andconfiguration of network topology and routing. Though these may revealtraffic anomalies or congestion problems, they do not in general revealpotential solutions. For instance, link load measurements may revealcongestion on a link, but shed little light on its cause, which ingeneral requires understanding the traffic matrix.

The inability of network operators to measure the traffic matrix is afundamental obstacle to developing sounds methods for network andtraffic engineering in operational IP networks.

SUMMARY OF INVENTION

The present invention is directed to mechanisms for measuring trafficvolume from a plurality of ingress points to a plurality of egresspoints in a large scale network, such as an IP backbone network. Thepresent invention infers the traffic matrix advantageously from widelyavailable link load measurements, such as SNMP data. First, such data iscollected and used to construct a gravity model of link to link trafficto capture an overall distribution of the volume of traffic. Additionalinformation on routing between points of ingress and egress for trafficflows can be incorporated to obtain significantly improved results,e.g., the model can incorporate information to model traffic exchangedwith peer networks in a typical IP backbone network. Second, the trafficmatrix is estimated by determining the matrix that minimizes a distancemetric to an initial tomographic solution based on the gravity model,subject to tomographic constraints. Quadratic programming can beutilized to determine the solution in the space of those admitted by thetomographic constraints closest to the solution obtained by the gravitymodel. This step advantageously does not require (higher-order)statistics or additional traffic modeling assumptions. Applying networkconfiguration and routing data to remove empty demands from the trafficmatrix serves to dramatically decrease the problem dimension ofcomputing the pseudo-inverse of the routing matrix. Then iterativeproportional fitting can be used to ensure that non-physical results arereplaced. The resulting traffic matrix comprises elements thataccurately estimate the traffic volume from the plurality of ingresspoints to the plurality of egress points in the network.

The present invention is especially accurate for large elements and isrobust, easily coping with data glitches and loss. The present inventionis fast and flexible, easily extending to incorporate more detailedmeasurements where available. The present invention enables true networkengineering, such as reliability analysis, traffic engineering, andcapacity planning, in an IP framework.

These and other advantages of the invention will be apparent to those ofordinary skill in the art by reference to the following detaileddescription and the accompanying drawings.

BRIEF DESCRIPTION OF DRAWINGS

FIG. 1 is an abstract diagram of an illustrative IP backbone network.

FIG. 2 is an abstract diagram of three points of presence in theillustrative IP backbone network.

FIG. 3 is a flowchart of processing performed to determine a trafficmatrix, in accordance with an preferred embodiment of an aspect of thepresent invention.

FIG. 4 is a flowchart of processing performed to remove the emptydemands in the initial tomographic solution, thereby reducing thecomplexity of the tomographic problem, in accordance with an embodimentof another aspect of the invention.

FIG. 5 is an graph abstractly illustrating how the gravity modelsolution is utilized to obtain the final solution using a tomographicapproach.

FIG. 6 is illustrative MATLAB source code for computing a weightedleast-squares estimate of the traffic matrix.

DETAILED DESCRIPTION

FIG. 1 is an abstract diagram of an illustrative Internet Protocol (IP)backbone network. The IP backbone network 100, as depicted in FIG. 1,can be represented as a set of nodes and links, associated with the IProuters 101, 102, 103, 104, etc., and the IP adjacencies between thoserouters. The nodes and links wholly internal to the network 100 arereferred to herein as backbone nodes and links; the others are referredto as edge nodes and links. A typical IP network managed by an InternetService Provider (ISP) will connect via edge links to other autonomoussystems (e.g., to a public Internet exchange point or directly to aprivate peer or transit provider) and customers (e.g., to a modem bankfor dial-up users, a web-hosting complex, or a particular business oruniversity campus). It is helpful to further categorize the edge linksinto access links 150, which connect to customers 160, and peering links140, which connect to other non-customer autonomous systems, e.g., peersA 110 and B 120, as depicted in FIG. 1. It is also helpful, without lossof generality, to characterize those routers that terminate accessand/or peering links as an edge router (ER) and those that terminateonly backbone links as a backbone router (BR). FIG. 2 furtherillustrates this terminology. Three points-of-presence (PoPs) 201, 202,203 in the IP backbone network are depicted, each with a number of BRs221, 222, 223, 224, 225, 226 and ERs 231, 232, 233, 234, 235, 236. Thelinks between BRs are referred to as core links, while the links betweena BR and an ER are referred to as non-core links.

It is assumed that there is access to routing information in the IPbackbone network 100. For example, in large IP networks, distributedrouting protocols are used to build the forwarding tables within eachrouter. It is possible to predict the results of these distributedcomputations from data gathered from router configuration files. Seeco-pending commonly-assigned U.S. Utility patent application Ser. No.09/876,384, entitled “TRAFFIC ENGINEERING SYSTEM AND METHOD”, which isincorporated by reference herein. There may be more than one routebetween two routers even using only shortest paths. It is assumed hereinthat traffic will be evenly distributed across all such routes, althoughthe present methodology can be easily adapted by one of ordinary skillin the art to handle unequal distributions.

It is also assumed that some form of link load measurements areavailable in the IP backbone network 100. Although the presenttechniques can work with different types of network load data, includingflow data, the present invention advantageously can take advantage ofdata widely available via the Simple Network Management Protocol (SNMP).See, e.g., D. Harrington, et al., “An Architecture for Describing SNMPManagement Frameworks,” Internet Engineering Task Force, Request forComments 2571, April 1999. Since every router maintains a cyclic counterof the number of bytes transmitted and received on each of itsinterfaces, it is possible to obtain basic traffic statistics for theentire network with little additional infrastructure support. SNMP datanotably has many limitations, which, preferably, should be taken intoaccount in implementing a useful methodology. Data may be lost intransit, either due to SNMP utilizing unreliable UDP transport or todata loss when copying to a central archive. Data may be incorrect,e.g., through poor router vendor implementations. The data collectionsampling interval is typically course, e.g. 5 minutes. Many of thetypical problems in SNMP data may be removed with minimal artifactsusing simple techniques. For instance, using hourly traffic averages,with five minute or finer data polls, mitigates the effect of missingdata substantially. Slightly more sophisticated methods of anomalydetection and interpolation can produce even better results. It isadvantageous to retain some of the unaveraged data for brief periods fortroubleshooting and alarming.

The traffic matrix can be computed with different levels of aggregationat the source and destination endpoints, for instance, at the level ofPoP to PoP, or router to router, or link to link. Intuitively, theaggregation level needs to be sufficiently high so that the trafficexchanged between different locations is not sensitive to the detailedcomposition of the traffic. On the other hand, when the aggregationlevel is too high (e.g., PoP to PoP), ISP routing policies operating ata more granular level may have a profound impact and can introduceserious systematic distortion to the overall traffic pattern. Theinventors have found it advantageous to utilize router to routermatrices, which are appropriate for a number of network and trafficengineering applications—and which, where useful, can be used toconstruct more highly aggregated traffic matrices (e.g., PoP to PoP)using routing information. Given two ERs E_(i) and E_(j), the trafficbetween these edge routers T_(ij) ^(E) is defined as the total amount oftraffic that enters the network at E_(i) and exits at E_(j), withT^(E)={T_(ij) ^(E)} the associated matrix. Similarly, given two BRs, thetraffic between these backbone routers T_(ij) ^(B) is defined such thatthe elements of the associated matrix, T^(B)={T_(ij) ^(B)}, refer totraffic that enters and leaves the core. The traffic matrix t will oftenbe referred to below in a vector form, in which the indices of thevector refer to source/destination pairs.

FIG. 3 is a flowchart of processing performed to determine the networktraffic matrix, in accordance with a preferred embodiment of the presentinvention. The processing can be readily performed by software orfirmware on a fast computer processor, such as a Sun Microsystems 336MHz UltraSPARC-II processor, with access to the above-mentioned linkload measurements and routing information.

At step 301, the link load measurements are collected. Where, forexample, SNMP data is being utilized, an SNMP “poller” can be readilyimplemented in the network to periodically request an appropriateManagement Information Base (MIB) from the relevant devices in thenetwork. The data can then be stored in an archive and averaged,aggregated, or interpolated, where appropriate, before furtherprocessing.

At step 302, the overall distribution of traffic in the network iscaptured using a simple model of the link-to-link traffic. For example,one of simplest approaches to modeling the link-to-link traffic is toapply what is referred to in the art as a “gravity model”. A generalformulation of a gravity model is given by the following equation:

$X_{i,j} = \frac{R_{i} \cdot A_{j}}{f_{i,j}}$where X_(i,j) is the matrix element representing the force from i to j;R_(i) represents the repulsive factors that are associated with“leaving” from i; A_(j) represents the attractive factors that areassociated with “going” to j; and f_(i,j) is a friction factor from i toj. Gravity models, which take their name from Newton's law ofgravitation, are commonly used by social scientists to model themovement of people, goods, or information between geographic areas. Inthe present context, X_(ij) can be interpreted as the traffic volumethat enters the network at location i and exits at location j; therepulsion factor Ri as the traffic volume entering the network atlocation i, and the attractivity factor A_(j) as the traffic volumeexiting at location j. The friction matrix, f_(i,j), can encode thelocality information specific to different source-destination pairs. Itis necessary to approximate the actual friction matrix using models withfewer parameters. In the embodiment described herein, a common constantis assumed for the friction factors, arguably the simplest among allpossible approximation schemes. The resulting network gravity modelsimply states that the traffic exchanged between locations isproportional to the volumes entering and exiting at those locations.

A simple gravity model can be utilized to estimate the amount of trafficbetween edge links as follows: Denote the edge links by l₁, l₂, . . . .The volume of traffic T(l_(i), l_(j)) that enters the network at edgelink l_(j) and exits at edge link l_(j) is estimated. Let T_(link) ^(in)(l_(i)) denote the total traffic that enters the network via edge linkl_(i), and T_(link) ^(out)(l_(i)) the corresponding quantity for trafficthat exits the network via edge link l_(i). The gravity model can thenbe computed by either of

${{T\left( {l_{i},l_{j}} \right)} = {{T_{link}^{in}\left( l_{i} \right)}\frac{T_{link}^{out}\left( l_{j} \right)}{\sum\limits_{k}\;{T_{link}^{out}\left( l_{k} \right)}}}},{{T\left( {l_{i},l_{j}} \right)} = {\frac{T_{link}^{in}\left( l_{j} \right)}{\sum\limits_{k}\;{T_{link}^{in}\left( l_{k} \right)}}{T_{link}^{out}\left( l_{j} \right)}}}$The first equation states that the traffic matrix elements T(l_(i),l_(j)) are the product of the traffic entering the network via edge linkl_(i) and the proportion of the total traffic leaving the network viaedge link l_(j), while the second is reversed and is identical undertraffic conservation—that is, the assumption that the interior networkis neither a source, nor sink of traffic. While this assumption isviolated (e.g., protocols running in the network interior act sourcesand sinks, and packet drops act as sinks) the volumes involved areinsignificant in the network considered. Most notably the actual resultsfrom the two equations are almost identical.

It is possible and preferable to generalize the above simple gravitymodel to take into account a wide range of additional informationprovided by link classification and routing policy. For example,typically, in large North-American ISPs, the majority of traffic isexchanged between network customers and network peers. The pattern andhandling of customer and peer traffic are quantitatively different, andthis has a large impact on the traffic matrix. Furthermore, this peeringtraffic has a large impact on every aspect of network design andengineering, and so estimating the associated traffic matrices is veryimportant. It is advantageous to adapt the gravity model to specificallydifferentiate between customer and peering traffic.

A generalized gravity model can be constructed as follows. It is assumedthat the network has a set of peers labeled P₁, P₂, . . . , andexchanges traffic with peer P_(i) over a set of edge links dedicated tothis peer. This is commonly termed private peering and is the dominantmode of peering for large IP backbones today. A set of customer accesslinks are labeled a₁, a₂, . . . , and a set of peer links are labeledp₁, p₂, . . . . The set of edge links carrying traffic to peer p_(i) isdenoted by

, and the set of all peer links by

. The set of all customer access links is denoted by

. SNMP measurements provide volumes of traffic on all edge links, j:T_(link) ^(in,out)(j), where the superscripts in (out) denotes trafficinto (out of) the backbone. The traffic entering, or exiting the networkto peer p_(i), is

${{T_{peer}^{x}\left( P_{i} \right)} = {\sum\limits_{p \in P_{i}}\;{T_{link}^{x}(p)}}},$where x=in or out. Then, the outbound traffic from access link α_(i)ε

to peering link p_(m)ε

_(j) is

${T_{outbound}\left( {a_{i},p_{m}} \right)} = \left\{ \begin{matrix}\frac{T_{link}^{in}\left( a_{i} \right)}{\sum\limits_{a_{k} \in A}\;{T_{link}^{in}\left( a_{k} \right)}} & {{T_{peer}^{out}\left( P_{j} \right)},} \\\; & {{{{if}\mspace{14mu} p_{m}} = {X\left( {a_{i},P_{j}} \right)}},} \\{0,} & {{otherwise}.}\end{matrix} \right.$where X is the exit peering link, as explained below. The inboundtraffic from peering link p_(i) to access link α_(j) is

${T_{inbound}\left( {p_{i},a_{j}} \right)} = {{T_{link}^{in}\left( p_{i} \right)}{\frac{T_{link}^{out}\left( a_{j} \right)}{\sum\limits_{a_{k} \in A}\;{T_{link}^{out}\left( a_{k} \right)}}.}}$

The internal traffic from access link α_(i) to access link α_(j) is

$\begin{matrix}{{{T_{internal}\left( {a_{i},a_{j}} \right)} = {\frac{T_{link}^{in}\left( a_{i} \right)}{\sum\limits_{{ak} \in A}\;{T_{link}^{in}\left( a_{k} \right)}}{T_{internal}^{out}\left( a_{j} \right)}}},{where}} \\{{T_{internal}\left( a_{j} \right)} = {{T_{link}^{out}\left( a_{j} \right)} - {\sum\limits_{p_{k} \in P}\;{T_{inbound}\left( {p_{k},a_{j}} \right)}}}} \\{= {{T_{link}^{out}\left( a_{j} \right)}{\left( {1 - \frac{\sum\limits_{p_{i} \in P}\;{T_{link}^{in}\left( p_{i} \right)}}{\sum\limits_{a_{k} \in A}\;{T_{link}^{out}\left( a_{k} \right)}}} \right).}}}\end{matrix}$These equations are developed from the following assumptions, whichreflect dominant Internet and ISP routing policies:

Transit peer (peering link to peering link) traffic. It is assumed thatthe volume of traffic that transits across the backbone from one peernetwork to another is negligible. Transit traffic between peers mayreflect a temporary step in network consolidation following an ISPmerger or acquisition, but should not occur under normal operatingcircumstances.

Outbound (access link to peering link) traffic. The proportionalityassumption underlying gravity modeling is applied on a peer-by-peerbasis: that is, the traffic exiting a specific peer comes from eachaccess link in proportion to the traffic on that access link. It isassumed that all of the traffic from a single access link to the givenpeer exits the network on the same peering link (as determined by theInterior Gateway Protocol (IGP) and Border Gateway Protocol (BGP)routing configuration). The exit peering link for traffic from accesslink α_(i) to peer P_(j) is denoted by X(α_(i), P_(j)), which may bederived from routing configuration information. The assumption istypically true in practice, except for example when short-term loadbalancing is performed. In such situations, the above methodology couldbe supplemented with available statistics on the affected prefixes,though the inventors' experience has been that the impact is small anddoes not affect the accuracy of the traffic matrix computation.

Inbound (peering link to access link) traffic. A network operator haslittle control over the injection of traffic into its network from peernetworks. Accordingly, it is assumed that the traffic entering from agiven peering link is split amongst the access links in proportion totheir outbound traffic.

Internal (access link to access link) traffic. It is assumed that thefraction of internal traffic from a given access link α_(i) to a secondaccess link α_(j) is proportional to the total traffic entering thenetwork at α_(i), and that the traffic between the links can be computedby normalization.

The generalized gravity model has been found by the inventors to matchactual Internet data very well. One possible explanation for this isthat geographic locality is not a major factor in today's Internet, ascompared to ISP routing policies. As long as the gravity model capturesthe essence of the routing policies, it becomes very accurate and thechoice of the above-mentioned friction factor is less critical. It iscertainly possible to further improve the above method by using a moreaccurate model with additional parameters. However, the margin forimprovement may be limited.

With reference again to FIG. 3, the result of the link-to-link gravitymodel at step 302 is a traffic matrix that is not guaranteed to beconsistent with the internal link measurements in the network—asignificant drawback. This is remedied, in accordance with an aspect ofthe invention, by utilizing a tomographic approach. It is not expectednor required that the above gravity model accurately model the trafficbetween all source-destination pairs. In fact, one would naturallyexpect certain pairs of locations to stand out from the overalldistribution, simply due to their specific characteristics (e.g., goingthrough transoceanic links). Rather, the model need only capture theoverall distribution. It is expected that tomographic estimation can beutilized to correct most of the violations in the assumptions underlyingthe model and thus significantly improve the accuracy.

Tomographic methods are based on the system of linear equations x=A twhere t is the traffic matrix, x represents the link loads, and A thenetwork routing matrix. The traffic matrix t is written as a columnvector t=(t₁, t₂, . . . , t_(m))^(T), where the M traffic matrixelements, t_(r), are the traffic between the rth source/destinationpair. The link traffic is the sum of the traffic matrix elements thatare routed across that link. The traffic (as measured in packets orbytes) that traverses the L links of the network during some period isrepresented as the set of observables x=(x₁, x₂, . . . , x_(L))^(T). TheL×M routing matrix A={A_(ir)} where

$A_{ir} = \left\{ {\begin{matrix}{F_{ir},} \\{0,}\end{matrix}\begin{matrix}{{if}\mspace{14mu}{traffic}\mspace{14mu}{for}\mspace{14mu} r\mspace{14mu}{traverses}\mspace{20mu}{link}\mspace{14mu} i} \\{otherwise}\end{matrix}} \right.$where F_(ir) is the fraction of traffic from source/destination pair rthat traverses link i. In essence, the equation x=A t states that thetraffic matrix must be consistent with network routing and measured linkloads throughout the network, not just at the edge. This matrix equalityis, however, highly under-constrained, and so allows many solutions.Tomographic methods differ in how a single “best” solution is identifiedfrom the possibilities. The majority of existing statistical tomographicapproaches (commonly referred to as “network tomography” methods) usemodels of the higher order statistics of the link load data to createadditional constraints. The present invention advantageously does notincorporate additional constraints, but rather uses the gravity model toobtain an initial estimate of the solution, which is further refined tosatisfy the constraints.

With reference again to FIG. 3, at step 303, the results of step 302 aretransformed into an initial traffic matrix solution, referred to hereinas t_(g). It is advantageous and preferable to transform thelink-to-link model results into a more tractable problem, such as aBR-to-BR traffic matrix using routing information. BR-to-BR trafficmatrices are generally more useful for traffic engineering tasks such asload balancing and are absolutely necessary for link/router failureanalysis. In addition, more highly aggregated traffic matrices, such asPOP-to-POP traffic matrices, may be directly derived from BR-to-BRmatrices using routing information. Even considering only a hundredbackbone routers, this leads to a problem with over a thousand unknowns,which is orders of magnitude more than the available constraints on linktraffic volume.

Accordingly, at step 304, the complexity of the tomographic problem forthe border router traffic matrix is advantageously reduced, in oneembodiment, by removing empty demands from the initial traffic matrixsolution. First, since there may be multiple BRs within each PoP,traffic will flow only between the closest of these as determined by IGProuting—thereby, rendering many of the BR-to-BR matrix elements empty.As depicted in the simplified illustrative topology of FIG. 2, there aretwo BRs in each PoP, connecting ERs within the PoP with redundant links.Given shortest path routing (and equal link weights on backbone links),it can be seen that all of the traffic from PoP B 202 to PoP C 203 willtraverse the route through BRs 222 and 223, while there will be notraffic entering the backbone nodes at BR 221 and departing at BR 224.While this is a very simple example, in operational IP networks, the setof paths consistent with IP routing will typically be significantly lessthan the set of all paths between router pairs.

FIG. 4 is a flowchart of processing performed to remove the emptydemands in the initial traffic matrix solution, in accordance with anembodiment of another aspect of the invention. The traffic matrix fromBR B_(i) to B_(j) is denoted T_(ij) ^(B). At step 401, all elements ofthe BR to BR traffic matrix are initially marked as empty. At step 402,the routing protocol, e.g., the Interior Gateway Protocol (IGP), issimulated to find the shortest paths between each source and destinationrouter. At step 403, for each path, let B_(i) and B_(j) be its first andlast BRs respectively, and mark T_(ij) ^(B) as not empty. Then, at step404, remove all T_(ij) ^(B) that remain empty. This step is equivalentto removing elements from t that will be zero because the correspondingroute is not used (unless failures occur).

Step 402 in FIG. 4 would typically be done for each pair of edgerouters, which can be prohibitive due to the large number of routers.However, the “topological equivalence” of edge routers can be exploitedto avoid having to run simulations for all possible pairs of routers.Two edge routers are said to be “topologically equivalent” if theyconnect to the same (non-empty) set of BRs and the protocol weights onthe corresponding links are the same. (It should be noted that there maybe multiple layer hierarchy of ERs within a PoP, and lower layer ERspass traffic to BRs only through the higher layer ERs. The lower layerof ERs may be removed from consideration as their (network compacting)traffic can be seen at the higher layer.) Such equivalent edge routerscan be grouped together and considered as what the inventors refer to asa single edge router equivalence class (EREC). The routes between thecomponent ERs of the same pair of ERECs advantageously should be thesame except for the first and last links. Consequently, only one IGPsimulation need be run for each pair of ERECs. Computing routes on theabove basis can reduce the computational burden by a factor of 20. Aftereliminating all of the empty demands, the number of unknowns can bereduced by a factor of 10, thereby turning a highly under-constrainedlinear inverse problem (with potentially hundreds of constraints andtens of thousands of unknowns) into a moderately under-constrainedproblem (with hundreds of constraints and about a thousand unknowns),and making the computation orders of magnitude faster.

Finally, at step 305 in FIG. 3, an optimization-based tomographicapproach is utilized to refine the initial gravity model solution. Asdepicted abstractly in FIG. 5, the final estimate of the traffic matrixis selected by choosing the solution in the space of those admitted bythe tomographic model that is “closest” to the solution obtained by thegravity model, in accordance with some form of objective function ordistance metric. For example, the gravity model may be refined using aleast-squares solution that minimizes the Euclidean distance to thegravity model solution subject to the tomographic constraints, asdepicted in FIG. 5. More specifically, the following quadratic programcould be solved

-   -   min ∥t−t_(g)∥    -   s.t ∥At−χ∥ is minimized        where ∥.∥ is the L₂ norm of the vector, i.e. the Euclidean        distance to the origin). As another example, a weighted        least-squares solution could be utilized, in which the        projection onto the subspace is not orthogonal, but rather        weighted by a function of the size of the estimated traffic        matrix elements (t_(g)). That is, the equation ∥(t−t_(g))/w∥        could be used as the objective function to minimize in the above        quadratic program, where w is the weight vector, and / is the        element-by-element vector division. As illustrated in FIG. 5,        the simple least-square solution is just an orthogonal        projection of the gravity model solution onto the constraint        sub-space. The weighted least squares solution, on the other        hand, gives different weight to different unknowns in the        solution.

Note that the tomographic constraints may be ill-posed due to possibledependency among different link constraints. Furthermore, theconstraints may not be satisfiable due to error and noise in the linkload data or possible routing changes that are not captured by thetopology data. The standard technique for dealing with ill-posedquadratic programs is to use Singular Value Decomposition (SVD) of therouting matrix A to compute its pseudo-inverse. The resulting solutionis the closest to the initial solution t_(g) among all solutions thatminimize the discrepancy against the tomographic constraints (∥At−x∥).Routines to compute the pseudo-inverse are available in many numericalcomputing packages, e.g. MATLAB. FIG. 6 lists MATLAB source code thatimplements such an approach.

It should be noted that, where an approach such as least-squares isutilized, the methodology may result in negative values, which areclearly without physical meaning. Accordingly, at step 306 in FIG. 3,these non-physical values should be corrected. One approach to avoidingthese values is to view the problem as a constrained optimizationproblem. However, a simple iterative procedure provides a fast andeffective alternative. Specifically, Iterative Proportional Fitting(IPF) can be utilized to ensure non-negativity. See J. Cao, D. Davis, S.V. Wiel, and B. Yu, “Time0 varying network tomography: router linkdata,” J. Amer. Statist. Assoc., Vol. 95, No. 452, pp. 1063-1075 (2000),which is incorporated by reference herein. For the initial estimate, thetraffic matrix estimated above can simply be used, with zero replacingthe negative elements of the matrix. IPF then proceeds by successivelyrefining the estimate using the above-disclosed method. The computationis simpler and the initial condition is not complex, since higher orderstatistics of the process are not modeled.

The above embodiment of the present invention is remarkably fast, takingas little as 5 seconds on a 336 MHz UltraSPARC-II processor to compute abackbone router to backbone router traffic matrix on a tier-1 IPnetwork.

The complexity of the generalized gravity model is O(N²) in the numberof edge links being considered, but the number of operations per term issmall. The worst case complexity of the above quadratic program islinear in the number of unknowns (elements in the traffic matrix), andquadratic in the number of constraints. In practice, however, thecomplexity of singular value decomposition methods is generally lessthan this. For instance, the SVD used in MATLAB are usually better thanthis complexity estimate would indicate. In reality, this computationcan run significantly under 2 seconds. Computation of the generalizedgravity model for a complete network on the order of 1000 routers cantake less than 3 seconds. Computing routes by taking advantage of edgerouter equivalence classes can reduce the computational burden by afactor of 20. The time taken in computing the routing matrix dominatesall other aspects of the above embodiment, taking two to three minutes.It should be noted, however, that this cost can often by amortized overmultiple traffic matrix computations because the routing matrix need berecomputed only when the network topology changes. The method used toreduce the problem size can be performed as part of computing therouting matrix with a computational cost that is a very small marginalcost on top of computing the routing matrix itself.

The inventors have tried the above embodiment utilizing straightleast-squares and a range of different weighting schemes, includingconstant weighting, linearly proportional to the terms in the gravitymodel traffic matrix, and/or proportional to the square root of thegravity model. The inventors have generally found that (a) the simplegravity model is better than the raw least-squares approach, (b) thegeneralized gravity model is better than the simple gravity model, and(c) the best results come from using the generalized gravity model withthe weighted least-squares method using square root weights (though theimprovement over using other weightings is small). Moreover, theinventors have found that, even where there are errors in the trafficmatrix, the results can still be used for a range of operational taskssuch as detecting traffic. changes. The overall approach is robust tomeasurement errors on the observables.

The foregoing Detailed Description is to be understood as being in everyrespect illustrative and exemplary, but not restrictive, and the scopeof the invention disclosed herein is not to be determined from theDetailed Description, but rather from the claims as interpretedaccording to the full breadth permitted by the patent laws. It is to beunderstood that the embodiments shown and described herein are onlyillustrative of the principles of the present invention and that variousmodifications may be implemented by those skilled in the art withoutdeparting from the scope and spirit of the invention. For example, thedetailed description describes an embodiment of the invention withparticular reference to IP backbone networks. However, the principles ofthe present invention could be readily extended to other networkarchitectures. Such an extension could be readily implemented by one ofordinary skill in the art given the above disclosure.

1. A method of measuring traffic volume from a plurality of ingress points to a plurality of egress points in a network, the method comprising: constructing a gravity model of link to link traffic utilizing link load measurements on links between the ingress points and the egress points to capture an overall distribution of the volume of traffic between the ingress points and the egress points; and finding a traffic matrix that minimizes a distance metric subject to tomographic constraints to an initial tomographic solution based on the gravity model, the traffic matrix further comprising elements that specify the traffic volume from the plurality of ingress points to the plurality of egress points in the network.
 2. The method of claim 1 wherein the initial tomographic solution is modified to remove empty demands based on simulations of routing in the network.
 3. The method of claim 2 wherein simulations are run only on routes that are not topologically equivalent.
 4. The method of claim 3 wherein the ingress points and the egress points are border routers in an Internet Protocol (IP) backbone network.
 5. The method of claim 4 wherein the gravity model differentiates between customer and peering traffic in the IP backbone network.
 6. A computer readable medium comprising executable instructions, the instructions when executed causing a processor of a computer to perform a method of measuring traffic volume from a plurality of ingress points to a plurality of egress points in a network, the method comprising: constructing a gravity model of link to link traffic utilizing link load measurements on links between the ingress points and the egress points to capture an overall distribution of the volume of traffic between the ingress points and the egress points; and finding a traffic matrix that minimizes a distance metric subject to tomographic constraints to an initial tomographic solution based on the gravity model, the traffic matrix further comprising elements that specify the traffic volume from the plurality of ingress points to the plurality of egress points in the network.
 7. The computer readable medium of claim 6 wherein the initial tomographic solution is modified to remove empty demands based on simulations of routing in the network.
 8. The computer readable medium of claim 7 wherein simulations are run only on routes that are not topologically equivalent.
 9. The computer readable medium of claim 8 wherein the ingress points and the egress points are border routers in an Internet Protocol (IP) backbone network.
 10. The computer readable medium of claim 9 wherein the gravity model differentiates between customer and peering traffic in the IP backbone network. 